'''
Created on Jan, 19 2011

@author: oabalbin
'''
import logging
import os
import sys
import subprocess
import tempfile

from optparse import OptionParser
from exome.jobs.base import JOB_SUCCESS, JOB_ERROR
from exome.jobs.job_runner import qsub_cac, qsub_loc, run_local
from exome.jobs.config import ExomePipelineConfig, ExomeAnalysisConfig, ExomeLaneConfig

# Global Variables
NODE_MEM=45000.0
NODE_CORES=12
SINGLE_CORE=1
MEM_PER_CORE= int(float(NODE_MEM) / NODE_CORES)
# wt=walltime
WT_SHORT= "24:00:00"
WT_LONG= "60:00:00" 


def read_length_trimmer(input_fastqc_file,output_fastqc_file, bases_to_trim, path_to_fastx):
    '''
    It uses the fastx toolkit to trim the reads 
    to a particular length. 
    bases_to_trim means the number of bases from the
    end of the read.
        [-h]         = This helpful help screen.
       [-f N]       = First base to keep. Default is 1 (=first base).
       [-l N]       = Last base to keep. Default is entire read.
       [-t N]       = Trim N nucleotides from the end of the read.
                      '-t'  can not be used with '-l' and '-f'.
       [-m MINLEN]  = With [-t], discard reads shorter than MINLEN.
       [-z]         = Compress output with GZIP.
       [-i INFILE]  = FASTA/Q input file. default is STDIN.
       [-o OUTFILE] = FASTA/Q output file. default is STDOUT.
    '''
    fastx_command=os.path.join(path_to_fastx,'fastx_trimmer')
    args=[fastx_command, '-t', str(bases_to_trim), '-i', input_fastqc_file, '-o', output_fastqc_file]
    
    args= [a.replace(',',';') for a in args]
    command1 = ",".join(args).replace(',',' ').replace(';',',')
    
    args2=['mv',input_fastqc_file,input_fastqc_file+'.orig']
    args3=['mv',output_fastqc_file,input_fastqc_file]
    command2= ",".join(args2).replace(',',' ').replace(';',',')
    command3= ",".join(args3).replace(',',' ').replace(';',',')
    
    command=command1+'\n'+command2+'\n'+command3

    return command
    
    
def read_quality_trimmer(input_fastqc_file,output_fastqc_file, bases_to_trim, path_to_fastx):
    '''
        [-h]         = This helpful help screen.
       [-q N]       = Minimum quality score to keep.
       [-p N]       = Minimum percent of bases that must have [-q] quality.
       [-z]         = Compress output with GZIP.
       [-i INFILE]  = FASTA/Q input file. default is STDIN.
       [-o OUTFILE] = FASTA/Q output file. default is STDOUT.
       [-v]         = Verbose - report number of sequences.
                      If [-o] is specified,  report will be printed to STDOUT.
                      If [-o] is not specified (and output goes to STDOUT),
                      report will be printed to STDERR.
    '''
    fastx_command=os.path.join(path_to_fastx,'fastq_quality_trimmer')
    args=[fastx_command, '-t', str(bases_to_trim), '-i', input_fastqc_file, '-o', output_fastqc_file]
    
    args= [a.replace(',',';') for a in args]
    command = ",".join(args).replace(',',' ').replace(';',',')
    
    return command

        
def read_dynamicTrimmer(input_fastqc_file, output_fastqc_file, probcutoff, path_to_DynamicTrim):
    '''
    -p or-probcutoff #     probability value (between 0 and 1) at which base-calling error is considered too high (default; P = 0.05) or
    -h or -phredcutoff #     Phred score (between 0 and 40) at which base-calling error is considered too high (default; Q = 13)
    -b or -bwa     use BWA trimming algorithm
    '''
    DT_command=os.path.join(path_to_DynamicTrim,'DynamicTrim.pl')
    args=[DT_command, input_fastqc_file, '-p', str(probcutoff), '-b']
    args= [a.replace(',',';') for a in args]
    command1 = ",".join(args).replace(',',' ').replace(';',',')
    # move the old mate file to
    args2=['mv',input_fastqc_file,input_fastqc_file+'.orig']
    args3=['mv',input_fastqc_file+'.trimmed'+input_fastqc_file]
    command2= ",".join(args2).replace(',',' ').replace(';',',')
    command3= ",".join(args3).replace(',',' ').replace(';',',')
    
    command=command1+'\n'+command2+'\n'+command3
    
    return command
    

def read_trimmer(analysis, config_file, num_processors, jobrunfunc):
    
    bases_to_trim=str(19) #config.base_to_trim
    path_to_fastx=config.path_to_fastx
    my_email=config.email_addresses
    
    for sample in analysis.samples:
        for lane in sample.lanes:
            
            lname = lane.name.split('_')
            jobn=lname[1][:5]+'_'+lname[2]
            
            lane_file = os.path.join(sample.output_dir, lane.name + ".xml")
            analysis.lane_dir
            job = ExomeLaneConfig()
            job.from_xml(lane_file, analysis.lane_dir)
            
            for mate,fastq_file in enumerate(job.fastq_files):
                print mate, fastq_file
                fastq_trimmed_file=fastq_file.replace('.fq','_trim.fq')
                
                args=read_length_trimmer(fastq_file,fastq_trimmed_file, 
                                    bases_to_trim, path_to_fastx)
                
                jobidtm = jobrunfunc(jobn+'_trim', args, SINGLE_CORE, cwd=None, walltime=WT_SHORT, 
                                     pmem=MEM_PER_CORE, deps=None, stdout=None, 
                                     email_addresses=my_email)

        
if __name__ == '__main__':
        
    optionparser = OptionParser("usage: %prog [options] ")
    optionparser.add_option("-r", "--config_file", dest="config_file",
                            help="file with run configuration")
    optionparser.add_option("-a", "--analysis_file", dest="analysis_file",
                            help="file with experiment configuration")
    optionparser.add_option("--paired_samples", dest="paired_samples", action="store_true", default=False,
                            help="paired samples snv calling")
    optionparser.add_option("--local", dest="local", action="store_true", default=False)
    optionparser.add_option("--cluster", dest="cluster", action="store_true", default=False)
    optionparser.add_option("-p", "--processes", type=int, dest="num_processors", default=1)

    (options, args) = optionparser.parse_args()    

    config = ExomePipelineConfig()
    config.from_xml(options.config_file)
    analysis = ExomeAnalysisConfig()
    analysis.from_xml(options.analysis_file, config.output_dir)


    if not (options.local ^ options.cluster):
        print "Must set either --local or --cluster to run job"
    if options.local:
        jobrunfunc = run_local
    elif options.cluster:
        jobrunfunc = qsub
        
    read_trimmer(analysis, config, options.num_processors, jobrunfunc)

